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Abstract 

A 2D, hexagonal in geometry, statistical model of fracture is proposed. The model is based on the 
drying fracture process of the bamboo Guadua angustifolia. A network of flexible cells are joined by 
brittle junctures of fixed Young moduli that break at a certain thresholds in tensile force. The system 
is solved by means of the Finite Element Method (FEM). The distribution of avalanche breakings 
exhibits a power law with exponent —2.93(9), in agreement with the random fuse model [BC06 . 
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1 Introduction 

Guadua angustifolia is an Andean Bamboo that 
is widely used both as a construction material, 
for structural frames and walls ( |TR09| ) . and as 
raw material for handicrafts, furniture and vessels, 
throughout Colombia ( |Jan 81 , AV93 ). Moreover, 
Guadua forests play a crucial ecological role, fos- 
tering native species, helping to regulate the water 
cycle, consuming CO2, and injecting water into the 
soil. Nowadays, it is having a revival as construction 
material for building houses, |LVV + 03] due to re- 
cent advances in research about bamboo structural 
properties and behavior }LEM09j as well as struc- 
tural laminated bamboo |XSC + 09j . 

Guadua structure is optimized to withhold 
stresses along its axis, with a high Young modu- 
lus in this direction, and to bend elastically for de- 
formations perpendicular to it, as all bamboos do. 
To achieve this performance, all fibers are oriented 
along the culm axis, more of them being nearer the 
outer radius than the inner one |NSWP06] , Thus, 
cracks appear mostly along the axis as a conse- 
quence of the drying process, due to the avalanche 
of parenchyma cells collapse. 

The statistical models of fracture (SMF) de- 
scribe macroscopical cracks by following the failure 
of individual elements ANZ06 . An SMF can be 



tracked by means of computer simulations, which 
can exhaustively explore the statistical properties 
of the model. For instance, fiber bundle mod- 
els are ID analytically solvable SMFs describing 
the degradation and failure of materials. ( Pci26, 
IDan45l ISor89l IKRHH061 IPHGedp . Furthermore, 
fracture models have also used to investigate the 
fractal nature of the cracks MP84. ChiOO , the size 
distribution on impact fragmentation KH96] and 
the hydraulic fracture as used on the oil indus- 
try |GTMC07j . All these models consist on indi- 
vidual elements interacting by some known force 
with their neighbors, where the breaking of a sin- 
gle clement redistributes the load on the surround- 
ing ones, causing an avalanche of breaks. Many 
of them exhibit a phase transition that shows it- 
self in a power-law distribution of these avalanches 
|KH96| IHKPK08) IKMHH03| IKRHH06) . On the 
other side, the finite element method has been suc- 
cessfully implemented in the field of computational 
failure mechanics to study discontinuities and fail- 
ure dB0lI lJirT)0llOHP+r)4] . 

In the present work we propose a 2D statistical 
model of fracture, inspired in the geometry and the 
mechanical processes of fracture of an Andean Bam- 
boo, Guadua angustifolia. This model is explained 
in section [2j and resembles the geometry of the dry 
foam models HP SA + 99j . but with key differences: 
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the rule for breaking a juncture here is not random 
but chosen from the internal forces on the system, 
the individual hexagonal cells do not share walls 
and there is no coalescence of neighbouring cells. 
As in other SMF models, the breaking of a single 
juncture redistributes the stress on the surrounding 
ones; but here it is done not by following a simple 
rule but by computing the new equilibrium state of 
the whole system. Actually, we use the Finite El- 
ement Method (FEM) to model the cell walls and 
junctures and compute such equilibrium state. All 
cells are similar in shape, but breaking thresholds 
are chosen at random. We explore the distribution 
of avalanches and found that this model exhibits a 
power-law behaviour (an indication of a phase tran- 
sition), as other SMF models do. 



2 Guadua 
Model 



structure 



and 



Bamboo culms have two main structures: flexible 
parenchymatous tissue and stiff fiber bundles sur- 
rounding the transport vessels (Figure [2]) . From a 
structural point of view, the former provide flexi- 
bility while the latter carry most of the weight of 
the plant. During industrial drying, a careful pro- 
cess has to be set up in order to reduce the oc- 
currence of macroscopical cracks that emerge from 
microcracks at the parenchymatous tissue between 
the fibers |MA06j . 

Our model, from now on the Cell Network Model 
of Fracture (CNMF), is a 2D version of the mechan- 
ical response of the parenchymatous tissue to the 
stresses produced by the final drying stage, when 
all intracellular water is gone. The CNMF is a sta- 
tistical model of fracture, following the breaking of 
individual elements. In the real bamboo, cells are 
joined by the so-called tricellular junctures (made of 
pectic polysaccharides and calcium |WB07[ IJMOOj ) , 
i. e. three individual polymer bonds connecting 
the cell walls by pairs. The model resembles this 
structure with uniform hexagonal cells, with six cell 
walls each, ordered in a honeycomb structure and 
connected at the corners by triangles of individual 
junctures. 

The model is solved by the Finite Element 
Method (FEM). Both cell walls and junctures are 
modeled as linear elastic elements, called Plane 
Frame Elements (PFE) by the FEM literature 
|Seg84[ |Bat961 |Hug00l . Each PFEs has six degrees 
of freedom: two translational ones and a rotational 
one at each of two nodes. They are characterized by 
its elastic modulus E, cross sectional area A, area 
moment /, and length L. The nodal displacements 
of an element spaning from node i to node j are 



given by ufj c = {ui, v i: 4>i, %, «j, 4>j} (Figure^. In 
element (or local) coordinates, the element stiffness 
matrix ( |Seg84| ) is given by: 
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The representation of this matrix in global coordi- 
nates is reached through the transformation 
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K( e )=(T( e )) T K; ( ^T( e ). The equilibrium condi- 
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(2) 



tions for the whole structure are expressed by 
K-u = f n , 



(3) 



where u, K and f n are the displacement vector, the 
stiffness matrix and the nodal force vector for the 
whole structure, respectively, obtained by assem- 
bling all individual elements Seg84 . This expres- 
sion allows us to compute the displacements from 
the forces, or viceversa. 

Disorder is introduced into the system by choos- 
ing at random the breaking thresholds of each junc- 
ture. The Young modulus is fixed at literature pub- 
lished values |YFZY07] .The initial geometry is reg- 
ular, so the length L and cross sectional area A of 
all elements are fixed. As parenchymatous cell walls 
arc ductile and tricellular junctures are brittle, only 
junctures are allowed to break if their strength sur- 
passes a global failure threshold. 

3 Evolution 

In the final drying stage, that is after the fiber sat- 
uration point when all intracellular water is gone, 
guadua wood is usually part of a structure, and 
the boundary conditions are fixed even by stress 
or strain. Let us impose that all boundary nodes 
are fixed in space, i. e. with displacement zero. 
The drying process is modeled as a tensile nodal 
force on all elements that is proportional to the hu- 
midity change Ah. This change increases until a 
first juncture breaks. Then, the forces redistribute 
among the resting elements, eventually causing an 
avalanche of sucessive breakings. The simulation 
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Figure 1: A) Bamboo culms [Lin] . B)A small Guadua angustifolia board showing the fibers embedded 
into the parenchymatous tisse and their distribution: more dense at the outer than at the inner side. C) 
3D Bamboo culm structure. Fiber bundles align along the axis of the culm ([Lic98 ). D) 2D photography 
of the Guadua parenchyma ( |LVV + 03] ). 



includes a linear elastic shrink step, a breaking cri- 
tcrium and a non-linear avalanche step. 

Linear elastic shrink step: Drying caused 
stresses are applied following a two step process. 
In the shrink step, each cell wall suffers shrinking 
nodal forces along its axis and acting at its ends, 

f ( n ] ={-/„, 0,0, /„, 0,0} (4) 

with fn—EAahAh. In analogy with the strain 
caused by temperature WO90 . ah is a constant 
defining how much force is caused by a unit humid- 
ity change Ah. So, increasing Ah is like advancing 
in the drying process. These forces ensembles the 
global vector f n . In the equilibrium step, the global 
displacements u are computed from Eq. ([3| . 

Breaking criterium From the global displace- 
ment vector u, the element displacement vectors 
vS e ^ are constructed. The actual total force f^ e ' 
acting on a single element due to the contractions 
and displacements of the whole structure is given 
by /< e )=K0u( e ). If the magnitude of the longitu- 
dinal component of this vector is larger than the 
breaking threshold, the element will break. By tak- 
ing Ah as control variable, the secant method gives 
us the minimum value for a first breaking, that is 
the drying point when the avalanche starts. At this 
point, Ah is fixed for the whole avalanche and that 
juncture is broken. 

Nonlinear avalanche step. The breaking of a sin- 
gle juncture has two consequences. First, both the 
global stiffness matrix and the global nodal forces 
have to be updated. This is achieved by reassem- 
bling them from all elements but the broken one, 
to obtain K new and /„ ncw . Second, the breaking 
causes a new global displacement, M ne w> which is 



computed from Eq. ([3]) . This new displacement may 
cause in turn a new juncture to break, which causes 
a new displacement, and so on, until no more junc- 
tures break and the avalanche stops. 

4 Results 

Our model consists of an arrangement of 25 cells 
(Figure |2j. The boundary nodes are fixed. Figure 
[4] shows the number of avalanches as a function of 
the avalanche size, for different values of the center 
of the distribution of juncture breaking thresholds. 
The upper line (+) corresponds to a threshold dis- 
tribution centered at I0~ 5 in units of the area of 
the cross section of the juncture, A, times the mean 
value of the elastic modulus, E. The lower line (o) 
corresponds to a value of 1 in the same units. We 
averaged over different realizations (this number is 
also written into the graph) . The behaviour of the 
curves suggests that there is a threshold value where 
the distribution of avalanches' sizes follows a power 
law. Our best fit occurs at a threshold of 0.35AE, 
with a power law exponent —2.93(9). This value co- 
incides within the error bars with the reported value 
o f 3.0 ([PHTirMlHH94] 'l for the random fuse model 
f |Her90p . This will suggesst that the different ge- 
ometry does not play a mayor role in the breaking 
statistics. 
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Figure 2: Upper left, the plane frame element spanning between the nodes i and j, arbitrarely oriented 
by an angle 0. Each node has two translational and one rotational degree of freedom. Lower left, two 
contiguous cells. Right, structure of the CNMF (not at scale). The hexagons represent the cells and the 
junctures are arranged into triangles. Hashing represents fixed boundary conditions. 









Figure 3: Model evolution. Left(5 steps), a gap appears due to the breaking of some junctures at the 
center of the sample. Middle, 7 timesteps, right 17 timesteps. 
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